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ABSTRACT 

A strategy is devised for a semi-coherent cross-correlation search for a young neutron 
star in the supernova remnant SNR 1987A, using science data from the Initial LIGO 
and/or Virgo detectors. An astrophysical model for the gravitational wave phase is 
introduced which describes the star's spin down in terms of its magnetic field strength 
B and ellipticity e, instead of its frequency derivatives. The model accurately tracks the 
gravitational wave phase from a rapidly decelerating neutron star under the restrictive 
but computationally unavoidable assumption of constant braking index, an issue which 
has hindered previous searches for such young objects. The theoretical sensitivity is 
calculated and compared to the indirect, age-based wave strain upper limit. The age- 
based limit lies above the detection threshold in the frequency band 75Hz < i/ < 
450 Hz. The semi-coherent phase metric is also calculated and used to estimate the 
optimal search template spacing for the search. The range of search parameters that 
can be covered given our computational resources ('~ 10^ templates) is also estimated. 
For Initial LIGO sensitivity, in the frequency band between 50 Hz and 500 Hz, in 
the absence of a detected signal, we should be able to set limits of _B > 10^^ G and 
e < 10"*. 

Key words: gravitational waves — pulsars: general — stars: neutron — supernovae: 
individual (SNR 1987A) 



1 INTRODUCTION 



The Laser Interferometer Gravitational Wave Observa- 
tory (LIGO) achieved i ts design sensitivity during its 
fifth science run [S5; (| Abbott et all 12009(^ 1]. Analysis 
of S5 data is progressing well, with new upper lim- 
its being plac ed on the st r ength of various classes of 
burst sources llAbbott et all l2009dl : lAbadie e t aT '2010a|; 



Abbott et all 2010a) 



2008 



lAbbott et all 



stocha stic backgrounds 
2009al). compact 



Giamp anig 

^ binary sources 

( Abbott et al.ll2009el: lAbadie et allgOlOdIdh a n d continuous- 
wave sources l|Abbott et al.l l2009j |bl. l2010bl : lAbadie et al] 
l2010bl ). In some cases, the LIGO limits on astrophysical 
parameters beat those inferred from electromagnetic astron- 
omy, e.g. the maximum ellip ticity and internal magnetic field 
strength of the Crab pulsar (| Abbott et al.ll2008l . [2010bl ). Re- 
cently, an S5 search was completed which placed upper lim- 
its on the amplitude of r-mode oscillations of the neutron 



star in the supernova remnant Cassiopeia A (|Abadie et al.l 
l2010bl ). 

Aspherical, isolated neutron stars constitute one 
promising class of co ntinuous-wave source candidates 
(Ostrikcr fc Gunnlll969l ). The origin of the semi-pe rmanent 
quad ru pole in these objec t s can b e thermoelastic jMelatod 
^Od; bsho mirskv et all l2000l: iNavvar fc OwenI l2006l: 



Haskell 2008) o r hyd ro magnetic (|Bonazzola fc GourgoulhorJ 
T99a iCu tlcr 2005: 'Haskell et al.' '2008*: iHaskelll 1200 



Akgun fc Wasscrman , 2008 ; Mastrano 2010). Thermoelastic 



deformations arise due to uneven electron capture rates in 
the neutron star crust. A persistent 5% temperature gra- 
dient at the base of the crust produces a mass quadrupole 



moment of 



lO-' 



g cm 



(^ 



10 ^; lUshomirskv et akl 



i200Qf ). Hydromagnetic deformations, on the other hand, are 
produced by large internal magnetic fields, and misaligned 
magnetic and spin axes. For example, a neutron star 
with spin frequency 300 Hz and internal toroidal field 
Bt > 3.4 X 10^^ G has an ellipticity e ~ 10"'' {{Bt)/W^^ C) 
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The deformation of an ideal fluid star with 
an arbitrary magnetic field distribution and a barotropic 
equation of state can be computed, with ellipticitie s as high 
as 10 ~^ predicted for some configurations (,ffaskell et alj 
l2008h . Additionally, some exotic neutron star models (e.g. 
solid str ange quark stars) allow for ellipticities as large 
as 10^^ lOwenl [2005h . Accreting neutron stars in binary 
systems can form magnetic mountains, with e ^ 10~^ 
(Melatos fc Pavnell200a IVigelius fc Melatosll2009al ). Defor- 
mations of all kinds relax viscoelastically and resistively 
over time, so that young neutron stars are expected to 
be generally stronger gravitational wave emitters. For 
example, thermoelastic deformations relax on the thermal 
conduction time-scale (~ 10* yr), afte r the temperature 
gradient in the crust has swit ched off (|Brown fc BildstenI 
1 19981 : IVigelius" fc Melatodlioiol ). Magnetic mountains relajc 
as the accreted matter diffuses throug h the magnetic field 
on th e Ohmic time-scale lO'^-lO** yr (jVigehus fc Melatod 
l2009bl ). 

A coherent search for 78 known radio pulsars was per- 
formed on S3 and S4 LIGO and GEO 600 data. Upper limits 
on the ellipticities of these pulsars were obtain ed, the small- 
est be ing e 7.1 x lO"'^ for PSR J2124-3358 jAbbott et al.l 
12007 j ). More recently, a coherent search for 116 known pul- 
sars was carried out using data from both the LIGO and 
Virgo detectors, placing an upper limit of e < 7.0 x 10~* for 
PSR J2124-3358. 

The youngest isolated neutron star accessible to LIGO 
probably resides in the supernova remnant SNR 1987A. 
The coincident detection of neutrino bursts from the su- 
pernova by detectors all over the world confirmed the core- 
coUa pse event, strongly indicating the formation of a neutron 
star jAgliettaet al.lTl987j ; iHirata et all 1 19871 : iBionta et all 
Il987l : iBahcall et ahiligSTh FI Constraints have been placed 
on the magnetic field strength, spin perio d, and other 
birth properties of the putative neutron star (|Michellll994l : 
iQgel man fc AlpaillioO^ : see Section [2] for details. However, 
searches for a pulsar in SNR 1987A have yielded no con- 
firmed sightings; upper limits on its luminosity h ave been 
placed in the radio, optical and X-ray b ands (Pcrc ival et al.l 
1 19951 : iBurrows et al] I2OO0I : iManchesteij |2007D. An uncon- 
firmed detection of a transitory 467.5 H z optical/infra-red 
pulsat ion in SNR 1987A was reported bv lMiddleditch et al.l 
i|2000l) . 

The likely existence of a young neutron star in SNR 
1987A makes it a good targ e t for gravitationa, l wave searches 
llPiran fc Nakamural 1 198^ : iNakamural Il989l '). A coherent 
matched filtering search was carried out in 2003 with the 
TAMA 300 detector, searching 1.2 x 10^ hours of data from 
its first science run over a 1-Hz band centered on 934.9 Hz, 
assuming a spin-down range of (2-3) x 10"^" Hz s~^. The 
search yielded an upper limit on the wave strain of 5 x 10~^^ 



^ An unconfirmed correlation was also reported between data 
taken by tlie Mont Blanc and Kamioka neutrino detectors 
and gravitational wave detectors in Maryland and Rome 
l|Amaldi et al.lll989h . Taken at face value, tiiese observations are 
consistent witli a weak neutrino pulsar operating briefly during 
the core-collapse eve nt. However, a ser i ous fl aw in the original 
analysis was found bv lDickson fc Schutj l l 19951) . whose reanalysis 
led them to conclude that the correlations were not physically 
significant. 



(jSoidaet al.ll2003l ). An earlier matched filtering search was 
conducted using 10^ hours of data taken in 1989 by the 
Garching prototype laser interferometer. The latter search 
was carried out over 4-Hz bands near 2 kHz and 4 kHz, did 
not include any spin-down parameters , and yielded an uppe r 
limit of 9 X 10~^^ on the wave strain (jNiebauer et al.lll993l ). 

There are two main types of continuous-wave LIGO 
searches: coherent and semi-coherent. The former demand 
phase coherence between the signal and search template over 
the entire time series. Although sensitive, they are restricted 
to small observation times and parameter ranges as they 
are computationally intensive. Semi-coherent searches break 
the full time series into many small chunks, analyse each 
chunk coherently, then sum the results in coherently, trad- 
ing of f sensitivity for computational load. ISantostasi et al.l 
(|2003l ) discussed the detectability of gravitational waves 
from SNR 1987A, estimating t hat a coherent search based 
on the lMiddleditch et all (|200(]| ) spin parameters requires 30 
days of integration time and at least 10^^ search templates 
covering just the frequency and its first derivative. In reality, 
the task is even more daunting, because such a young object 
spins down so rapidly, that five or six higher-order frequency 
derivatives must be searched in order to accurately track 
the gravitational wave phase. A Bayesian Markov Chain 
Monte Carlo method was proposed as an alternative to 
cover the parameter space efficiently (' Umstatter et al]|2004l : 
Umstaetter et al. 2008 ). As yet, though, SNR 1987A has not 
been considered a feasible search target, because even Monte 
Carlo methods are too arduous. In this paper, we show how 
to reduce the search space dramatically by assuming an as- 
trophysically motivated phase model. 

In this paper, we discuss how to use a cross-correlation 
algorithm to search for periodic gravitational waves from a 
neutron star in SNR 19 87A. The search is semi-coherent 
(jPhurandhar et al] I2OO8I ) . The signal-to-noise ratio is en- 
hanced by cross-correlating two data sets separated by an 
adjustable time lag, or two simultaneous data sets from dif- 
ferent interferometers, thereby nullifying short-term timing 
noise (e.g. from rotational glitches). This is a modification 
of the metho d used in searches for acosmological stochastic 
background (jAbbott et al.l 2007al. l2009al) and f or the low- 
mass X-ray binary Sco X-1 ( Abbott et al.ll2007bl ). In Section 
[21 we review the properties of SNR 1987A and its putative 
neutron star. Section[3]briefly describes the cross-correlation 
algorithm and the data format. We estimate the theoretical 
sensitivity of the search in Section ID Section [5] describes 
an astrophysical model, which expresses the gravitational 
wave phase in terms of the initial spin, ellipticity, magnetic 
field, and electromagnetic braking index of the neutron star. 
We calculate the semi-coherent phase metric and the num- 
ber of templates required for the search in the context of 
the astrophysical phase model. Given the computational re- 
sources available to us, we derive upper limits on the grav- 
itational wave strain, ellipticity and magnetic field which 
can be placed on a neutron star in SNR 1987A with a cross- 
correlation search. Finally, SectionQsummarises the results. 



2 A YOUNG NEUTRON STAR IN SNR 1987A 

SNR 1987A is the remnant of a Type II core-collapse su- 
pernova which occurred in February 1987, 51.4 kpc away 
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in the Large Magellanic Cloud (a = Sh 35m 28 .03s, 5 
= -69° 16' 11.79" ) (see reviews by iPanagial boosi and in 
llmmlcr ct al."2007'). It s progenitor was the blue supcrgiant 
Sk 1 ( Panagia 1987: Gilmozzi et al1 [l987: Barkat & Whcclcr 
1 19881 : IWooslev et all l2002h . The color of the progeni- 
tor, as well as the origin of the complex three-ring 
nebula in the remnant, are still unexplained. Detailed 
simulations of the evolutionar y hist ory of Sk 1, per- 
formed by iPodsiadlowski et al.l (|2007l ). support the the- 
ory that two massive stars merged to form an oversized 
2OM0 red supergiant 2 x 10^ years before the supernova, 
which eventually shrank as its envelope evaporated (e.g. 
IPodsiadlowski fc Joss|[l989l : IPodsiadlowski et al.lll99d ). An 
alternative theory suggests that Sk 1 was instead a single 18- 
20 Mq red supergiant which evol ved into a blue supergiant 
via wind-driven mass loss (e.g. IWooslevI Il988l : i aio et al.l 
[l988. : Sugcr man et al.l liboS'). 

There is strong evidence for the existence of a neu- 
tron star in SNR 1987A. The progenitor mass range re- 
quired to produce Type II supernovae, 10-25 Mq , which in- 
cludes the above evolutionary scenarios, is the sa me range 
requi r ed to produce neu tron star remnants ( Wooslev et al.l 
l2002l : iHeger et al.1 12003|) . The secure neutrino detections 
mentioned in Section [1] support this conclusion. Although 
there have been no confirmed pulsar detections, numerous 
searches have placed upper limits on th e fiux and luminos - 
ity at radio (< 115 /iJy at 1390 M Hz. lManchesteij | 2007[) . 



optical/near-UV (< 8 x 10'^'' ergs s ^ , iGraves fc et al. 
and soft X-r ay (< 2.3 x 10^'* erg s~^, iBurrows et al 



200£), 



wavelengths. iMiddleditch et al] (2000) reported findini 



I2OOOI ) 



ig an 

optical pulsar in SNR 1987A with a frequency of 467.5 Hz, 
modulated sinusoidally with a ~ 1-ks period, consistent 
with precession for an ellipticity of e ~ 10~^. However, 
the p ulsations wore reported to have disappeared after 1996 
iMld dlcditch ct al. 2000) and were never confirmed indepen- 
dently. 

There are several possible reasons why a pulsar in SNR 
1987A has not yet been detected. If its spin period is greater 
than 0.1 s, it woul d not be b right enough to be detectable in 
the optical band (jPacini fc S alvati 1987; Manchest er 2007). 
If the radio emission is incoherent or the emission region is 
patchy, the pulses may have been missed, even if the beam 
width is as wide as is typical for yo ung pulsars l|Manchesterl 
l2007h . IShternin fc YakovlevI (|2008l ) argued that, although 
the neutron star's theoretical X-ray luminosity exceeds the 
observational upper limits by a factor of 20-100, the current 
upper limits still allow fo r concealment behind an opaque 
shell formed by fall back (IWoosley fc Weaver] IT995I). How- 
ever, simulations by iFrver et al. I l|l999l) suggest that, once 
fallback ceases, the accreted material cools, leaving no ob- 
scuring atmosphere. 

Another possible reason why a pulsar has not yet been 
detected is that its magnetic field is too weak. The weak- 
field theory is supported by theoretical models, in which the 
field grows only after the neutron st ar is formed and can 
take up to 10^ years to develop (e.g. iBlandford fc Romanil 
1 19881 : lReiseneggeill2003l). A grow th model for SNR 1987A 
was proposed by iMichell l| 19941 ). in which the magnetic 
field of a millisecond pulsar intensifies from 10^" G at 
birth to ~ 10^^ G after several hundred years (exponen- 
tial and linear growth were considered, yielding growth 
times of ~ 0.3-0.7 kyr), before the pulsar has time to spin 



down significantly. In an alternative model, the neutron 
star is born with a strong magnetic field, which is ampli- 
fied duri ng the first few seco nds of it s life by d ynamo ac- 
tion (e.g. iDuncan fc Thompson 1992; Bonanno et a l. 200^). 
Assuming this model, measurements of the known spin 
periods of isolated radio pulsars imply a distribution of 
birth magnetic field strengths between 10^^ G and 10^'^ G 



(|Arzoumanian et al.ll20o3 : Faucher-Giguere fc Kaspi '200f|). 



Several birth scenarios for the pulsar in SNR 1987A were 
considered bv lOgelman fc Alpaij (|2004l ) in this context, who 
concluded that the maximum magnetic dipole moment is 
< 1.1 X lO^'' G cm^ 2.5 x 10^* G cm^ and 2.5 x 10^° G cm^ 
for birth periods of 2 ms, 30 ms, and 0.3 s respectively. How- 
ever, the dynamo model also accommodates a magnetar in 
SNR 1987A, with magnetic dipole mom ent > 2.4 x 10^"^ G 
cm'^, regardless of the initial spin period (lOgelman fc Alp"arl 
I2OO4). 

Estimates of the birth spin of the pulsar in SNR 
1987A are more uncertain. Simulations of the bounce and 
post-bounce ph ases of core collapse were performed by 
lOtt et all (|2006t ) to determine the correlation between pro- 
genitor properties and birth spin. These authors found 
proto-neutron star spin periods of between 4.7-140 ms, 
proportional to the progenitor's spin period. A Monte 
Carlo popula tion synthesis study usin g known velocity 
distributions (jArzoumanian et al.l 12002 favoured shorter 
millisecond periods, but a sim ilar population study by 
iFaucher-Giguere fc Kaspil (|2006l ) argued that the birth spin 
periods could be as high as several hundred milliseconds. 
Faint, non-pulsed X-ray emission from SNR 1987A was first 
observed four months after the supernova and dec reased 
steadily in 1989 (|Dotani et al.ll 19871 : llnoue et al.lll991^ ■ lead- 
ing to the suggestion that a neutron star could be powering 
a plerion that i s partially obscured by a fragmented super- 
nova envelope. iBandiera et al.l (|l988l ) modelled the X-ray 
spectrum from a nebula containing a central pulsar, with a 
magnetic field of 10^^ G and an expansion rate of 5 x 10* cm 
s"^ The authors found a fit to the SNR 1987A data for a 
pulsar spin period of 18 ms. 



3 THE CROSS-CORRELATION ALGORITHM 

In this section, we bri efiy summarise the cross- correlation 
method described in iDhurandhar et al] (|2008l ). a semi- 
coherent search algorithm designed specifically to search 
for continuous-wave gravitational radiation. It operates on 
Short Fourier Transforms (SFTs) of data segments of length 
AT = 30 min, whose duration is chosen to minimise the 
Doppler effects due to Earth's rotation. In each SFT, the 
fcth frequency bin corresponds to the frequency V}^ — k/ AT 
for < fc s; N/2 and Uk = (fc-iV)/AT for N/2 fc < iV-1, 
where A'^ is the total number of frequency bins in the SFT. 

The output x{t) of a detector is the sum of the instan- 
taneous noise, n{t), and the gravitational wave signal, h{t). 
The noise is assumed to be zero mean, stationary, and Gaus- 
sian. Its power is characterised by Sn{iy), the single-sided 
power spectral density (i.e. the frequency-dependent noise 
fioor) in the following way: 



(1) 



where * denotes complex conjugation. Therefore, in the low 
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signal limit {\h(t)\ <^ \n{t)\), the power in the k-th frequency 
bin of SFT I can be approximated by 



(2) 



where we apply the finite time approximation to the delta 
function in ((TJ, i.e. 5at{v) = sm{T:vAT)/{nv) Ri AT. 

In the cross-correlation algorithm, SFTs are paired ac- 
cording to some criterion (e.g. time lag or interferome- 
ter combination) and multiplied to form the raw cross- 
correlation variable 



(3) 



(Ar)2 ' 

where I and J index the SFTs in the pair. The gravitational 
wave signal is assumed to be concentrated in a single fre- 
quency bin in each SFT (because AT <^v jv due to sidereal 
or intrinsic effects), whose index is denoted by ki or kj. 
The frequency bins in the two SFTs are not necessarily the 
same; they are related by the time lag between the pair and 
between interferometers, as well as spin-down and Doppler 
effects. For an isolated source, the instantaneous frequency 
at time t is given by 

v(t) = u(t) + v{t)^, (4) 
c 

where v{t') is the instantaneous frequency in the rest frame of 
the source, v is the detector velocity relative to the source, 
n is the position vector pointing from the detector to the 
source, and c is the speed of light. The instantaneous signal 
frequencies in SFTs / and J, vi and vj, are calculated at 
the times corresponding to the midpoints of the SFTs, Ti 
and Tj. The frequency bin kj is therefore shifted from kj by 
an amount ATSvij, with 5vij — uj — vi (|Dhurandhar et ahl 
[2008). For convenience, we now drop the subscripts ki and 
kj. 

In the low signal limit, yij is a random, complex vari- 
able. The cross-correlation statistic comprises a weighted 
sum of yij over all pairs IJ. yij has variance af j = 
Si'''(iyi)Sl;'\uj)/{4:AT^), where Si'\vi) is the power spec- 
tral density of SFT I at frequency vi, and SiJ\i'j) is the 
power spectral density of SFT J at frequency pj. 

The parameters describing the amplitude and the phase 
of the signal are contained within the signal cross-correlation 
function Gij, defined as 

+Fi^F.,^Al ^ i{Fi+F.j^ - Fiy,F.,+ )A+Ay] ,(5) 

with A^ij = ^i{Ti) - ^j{Tj). $i{Ti) and ui{Ti) are the 
phase and frequency at time Ti, whereas ^j{Tj) and vj{T,j) 
are evaluated at time Tj. Note t hat there is an error in equa- 
tion (3.10) of lDhurandhar et al.l (i200.S)i which omits the fac- 
tor of e-^'^'^'ri''^(^^>-''-'(^-''l arising from the choice of time 
origin of the Fourier transforms. The phase factors are deter- 
mined by the astrophysical phase model described in Section 

m 

The terms in square brackets in ((5]) depend on the po- 
larization angle ip, and the inclination angle t between n and 
the rotation axis of the pulsar, in the following way: 



A+ = 



1 -I- cos i 



COS t. 



(6) 
(7) 



_F + (t;n, ^) = a(t; n) cos 2'(/) -|- 6(f; n) sin 2-i/), 
Fx (t; n, ^) — fe(i; n) cos2'!/) — a(t; n) sin2'i/'. 



(8) 
(9) 



where a{t\ n) and b(t\ n) are the detector response functions 
for a given sky position, and are de fined in equations (12) 
and (13) of I Jaranowski et al.l lll998|). A geometrical defini- 
tion is also given in IPrix fc Whelanl(|2007l ). The gravitational 
wave strain tensor is 



h (t) = hoA+ cos $(t) V+ 4- hoA^ sin $(t) V 



(10) 



where ho is the gravitational wave strain, and tf +,x are the 
basis tensors for the + and x polarizations in the transverse- 
traceless gauge. 

In principle, one should search over the unknowns cos l 
and xp, but this adds to the already sizeable computational 
burden. Accordingly, it is customary to average over cos l 
and i/) when computing Qij, with 

{gij)..s.,^ = ^exp-^*^.' ^--^^['^^(^^'-''■^(^■^'l(a,a.+6,6.), 

(11) 

where ai^j — a{Ti^,j; n) and bi^j — b{Ti^,j\ n). Once the first- 
pass search is complete, a follow-up search on any promising 
candidates can then be performed, which searches explicitly 
over cos l and Preliminary Monte Carlo tests indicate that 
the detection statistic resulting from pip is approximately 
10—15% smaller than if the exact cos l and V" values are used. 

The cross-correlation detection statistic is a weighted 
sum of yij over SFT pairs. The number of pairs which can 
be summed over are limited by the available computational 
power. We discuss the computational costs of running the 
search in Section [6.21 The cross-correlation detection statis- 
tic is given by 



p = T.ij{uijyij 4- u'ljy'ij), 
where the weights are defined by 



UIJ 



(12) 



(13) 



For each frequency and sky position that is searched, we 
obtain one real value of p, which is a sum of the Fourier 
power from all the pairs. Ignoring self-correlations (i.e. no 
SFT is paired with itself), the mean of p is given by pp — 
^0 ^ij\Gij\'^ /'^'ij- In ths low signal limit, the variance of p 
is CTp = 2E/j|tJ/j|^/crf J. In the presence of a strong signal, 
and if self-correlations a re included, pp and cTp scale as 
(jPhurandhar et al.ll2008l ). 



4 SENSITIVITY 

4.1 Detection threshold 

Detection candidates are selected if they exceed a threshold 
value, pth - For a given false alarm rate Fa , this threshold is 
given by (jPhurandhar et al.|[200l ) 



pth = 2''''aperfc-'(2F,/iV), 



(14) 



where erfc is the complementary error function, and A'^ is 
the number of search templates used. In the presence of a 
signal, the detection rate for events with p > pth is given by 



^ ^ ierfc ( ^HL^ 



(15) 
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As Hp oc /ig, one can calculate the lowest gravita- 
tional wave strain that is detectable by the search to be 
ijPhurandhar et alJliooj ) 



hth{v) = 



5I/2 



1/4 

airs 



-,1/2 



AT 



(16) 



In we define S = erfc-^(2Fa)+erfc-^(2Fd), Fd is the 

false dismissal rate, is the mean-square of the signal 

cross-correlation function defined in ([5]), A^pairs is the number 
of SFT pairs, and Sn{v) is the single-sided power spectral 
density of the interferometers (assumed to be identical). 

One can estimate (IGul^)^^* theoretically for the spe- 
cial case where Tj — Tj and Qij is averaged over cost, ip, 
and sidereal time. In this case, the primary contribution to 
A<I>/j is the term [r(T/) — r(Tj)] • n/c, where r{t) is the 
position of the detector at time t in the frame of the solar 
system barycentre. Under these assumptions, equation pi|) 
can be expressed in terms of the overlap reduction func- 
tion l|WhelaDll2006l ). which depends only on v,a, and 5. For 
SNR 1987A, we have (a,S) = (1.46375 rad, -1.20899 rad), 
and hence (je/./i^)"^'''' = 4.6882. Assuming Fa = Fa = 0.1, 



AT = 1800 s, and A^pairs = 
SFTs), equation (fT6)) gives 



lO'' (approximately 1 year of 



hti,{u) = 5.92 X 10" 



1 1/2 



Hz- 



(17) 



Figure [T] is a graph of ^th as a function of i/. The val- 
ues of Sn{i^) axe based on LIGO's S5 noise characteristicsl3 
The S5 run began in November 2005 and accumulated a 
year's worth of triple coincidence data. For a signal from 
SNR 1987A to be detectable, we must have hth ^ ho- 



4.2 Minimum ellipticity and indirect, age-based 
limit 

The deformation of a neutron star is parameterised by its 
ellipticity e. The gravitational wave strain at Earth emitted 
by a biaxial neutron star is 



ho = 



(18) 



where G is Newton's gravitational constant, c is the speed 
of light, I is the moment of inertia, D is the distance to the 
source, and v is the gravitatio nal wave frequency, assu med 
to be twice the spin frequency jjaranowski et anil998h . 

An upper limit on ho can be derived from existing elec- 
tromagnetic data by assuming all the observed spin down 
comes from the gravitational wave torque, i.e. the observed 
frequency de rivative v satisfies v = —{32n'^Ge^ lu^) / {5c 



1038kg mVs 



1/2 



1/2 



D 



51.4kpc 



)(20) 



Unfortunately, without having observed any pulsations 
from SNR 1987A, it is impossible to determine f or |/>| a 
priori. Instead, we note that i> can be re-expressed in terms of 
the characteristic age of the source, Tc — !^/(4|/>|), assuming 
that 1/ today is much less than u at birth. The factor 4 
arises if one assumes the gravitational radiation dominates 
electromagnetic spin down, in order to remain consistent 
with (|19|) : in reality, electromagnetic spin down is expected 
to dominate, with Tc = i^/(2|!>|). Equation 1)20^ then reduces 
to 

-1/2 



hth{iy) ^ 3.39 X 10" 



D 



(21) 



_19yry V51.4kpc^ 

The right-hand side of (|21|) is graphed as a horizontal red 
line in Figure [T] The detectability condition (|2ip is then 
satisfied for spins in the range 75 Hz < J' < 450 Hz. Note 
that we have chosen Tc — 19 yr, the age of SNR 1987A in 
2006 when the S5 search began. 

It is important to note here that the assumption that f 
is currently much less than at birth is likely untrue for the 
object in SNR 1987A, as it is so young. Hence, the indirect, 
age-based limit in equation (I21|l and the horizontal line in 
Figure [T] are only indicative of the expected gravitational- 
wave emission strength (in fact, they are upper limits). Ex- 
act calculations of u and are performed in Section [5. II 



5 AN ASTROPHYSICAL MODEL FOR THE 
GRAVITATIONAL WAVE PHASE 

All continuous wave searches to date have used the stan- 
dard model for the gravitational wave phase, described in 
terms of a Taylor expansion i nvolving spin frequency deriva- 
tives (jjaranowski et al]|l998l ). For a young object like SNR 
1987A, which spins down rapidly, it is not computationally 
feasible to search over the six or more frequency derivatives 
typically needed to track the phase accurately. In this sec- 
tion, we present an alternative model for the gravitational 
wave phase, stated in terms of astrophysical parameters (i.e. 
the magnetic field strength and the neutron star ellipticity) 
instead of spin frequency derivatives. It tracks the phase ex- 
actly using four parameters, under the restrictive assump- 
tion (justified further below) that the braking index is con- 
stant. 

The phase of a slowly evolving gravitational wave signal. 



$(*) = 'l'(to) + 27r / dt,y{t), 

Jto 

(jWette et aLllioOa ). Combining this with mi to eliminate e can be approximat ed by the 



ho < 



1 



/ 5GI\i 



1/2 



(19) 



D \ 2c^v 

Hence, for SNR 1987A to be detectable (i.e. /ith ^ ho), we 
require 

/ith(i^) s; 1.66 X 10"^° 



l|jaranowski et aLlll998l ) 

s 

ii{t) = ^to) + 2'K^J' 



(22) 



Taylor expansion 



j-fc+i 



(fc + 1): 



27rn • vitn) -c-^ 

fe=0 



k\ 
(23) 

where t''*' is the fc-th derivative of the gravitational wave 
frequency at time to, and a is the number of spin-down pa- 
rameters required to achieve a given accuracy. The compu- 
tational cost of using (|23p is substantial for rapidly deceler- 
ating objects. For a maximum allowable phase error of one 



^ Available at'http://www.ligo.caltech.edu/jzweizig/distribution/LSC_B^tel(!, the maximum bin size in the fc-th derivative is ' is 
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Figure 1. Theoretical sensitivity of tlie cross-correlation search for SNR 1987A as a function of gravitational wave frequency (blue curve), 
assuming the initial LIGO detector power spectral density, a false alarm rate of 0.1, and a false dismissal rate of 0.1. The blue curve 
shows the theoretical sensitivity for the special case where the search uses 10^ pairs of time-coincident 30-minute SFTs, and averages 
over inclination angle, polarization angle, and sidereal time (see discussion in Section 14.111 . The horizontal red line shows the indirect, 
age-based limit assuming v <C i^b (see discussion in Section l4.2l l. 



Ai/('=' = (fc + l)!/Titg\ implying TVfe ^ i/C^'/Ai^e-'' templates 
in that derivative and A^totai = Hfe^o templates overall. 
We discuss this matter further in Section [6.21 

To improve on the above situation, we recognize that !> 
for an isolated neutron star is the sum of gravitational-wave 
and electromagnetic torque contributions: 



Q/ 5 ^/ n 
if - Q2I' , 



(24) 
(25) 



where is the neutron star radius, B is the polar magnetic 
field, n is the electromagnetic brak ing index (theoretically 
equa l to 3, but could be as low as 1 .8: iMelato j| 19971 : iPalonibal 
|2005| ). Assuming that the electromagnetic torque is propor- 
tional to a power of u, then 1/ must enter the torque in the 
combination Ri,v/c, (i.e. the ratio of ii* to the characteristic 
lever arm, the hght cylinder distance, c/2nv) on dimensional 
grounds. In terms of an arbitrary reference frequency, v,:ci, 
we write v = ~Qi (y jv^^if — Qi {u/vref)", with Qi = Qiv^^t 
and Q2 = Q'2^]^ci- Throughout this paper, we set Vy-d = IHz 
for simplicity. 

There may, of course, be other torques acting on a newly 
born neutron star. For example, nonlinear r-mode instabil- 
ities can emit a significant a mount of gravitat ional radia- 
tion under certain conditions (|Owen et al.lll99^ '). If there is 
a rapidly rotating pulsar with B 4, lO" G in SNR 1987A, 



its instability time scale (27 years) would exceed its age, 
and the gravitational radiation from the instabilities alone 
should be detectable by Advanced LIGO (Brink et al . 20ol: 
iBondarescu et al.l 12003 ). However, for the purposes of our 
search, we assume that the spin down is described by (|25p . 
An equally serious issue is that n may change over the 1 yr 
integration period, although in (|25|) . we assume that n is 
constant. Young pulsars have n < 3, and it can b e argued 
that n approaches 3 over the spin-down time-scale l|Melatosl 
In this search, we maintain the assumption of con- 
stant n. However, it is possible to extend (|25|) to include 
time-dependent n in future searches. We aim, in the first 
instance, to exclude the simplest astrophysical model while 
recognizing that it covers only a small fraction of the total 
parameter space. 

When implementing the search, instead of stepping 
through a grid of frequency derivatives, we search instead 
over v,Q\,Q2, and n. This reduces the number of parame- 
ters and allows one to track the phase more accurately for 
a given computational cost, as errors stemming from incor- 
rect choices of {v, Qi, Q2,n) grow more slowly with observa- 
tion time than errors stemming from higher-order frequency 
derivatives. The improvement is quantified in Section 16.21 
We note that the search targets a source with a known po- 
sition, hence in our estimates we consider only a single sky 
position. 
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5.1 Historical spin down 

We can use the possible spin histories of a source like SNR 
1987A with a known age to constrain the invisible values 
of (i^, Qi, Q2, n) today and hence the maximum amount of 
phase evolution to be expected during a LIGO integration. 

There are two ways of estimating p and u for a source 
whose age is known. In the simplest situation, where the 
current spin frequency v is much smaller than the birth 
frequency Vb, the characteristic age Tc ~ —v/ [{{n) — 1) z>] 
closely approximates the true age irrespective of Ub, where 
(n) is the mean braking index, averaged over the time since 
birth. Under these conditions, a source with unknown v and 
v lies on a line of slope — Tc ((n) — 1) in the v-v plane. How- 
ever, as discussed in Section [4.21 this is not necessarily true 
for SNR 1987A, which was only 19 years old at the start of 
the S5 search. In order to calculate v and v exactly with- 
out using the characteristic age approximation, one must 
integrate (|25|l over the lifetime of the source. Accordingly, 
we adopt this approach and map out the regions in the v-v 
plane which can be reached from Vb by electromagnetic-plus- 
gravitational- wave spin down and physically sensible choices 
of e, B and n. 

Figure [2] shows the range of possible u and z> values at 
f = 19yr obtained by solving for 10"" < e ^ 10"^, 
10"-^ s; B ^ 10^^ G and 0.1 kHz < i^b 1.2 kHz. For refer- 
ence, we plot the search sensitivity (black curve in the v-v 
plane) obtained from (|17|) . According to (|17[) . the search is 
only sensitive to combinations of v and u above the black 
line. The conservative limits set by the characteristic age ap- 
proximation are plotted as cyan lines. The lines correspond 
to (n) = 1.8 (top), (n) = 3 (middle) and (n) = 5 (bottom). 
For a given value of (n), an object lies on the line for u <^ Ub, 
and below the lines for v "^Vb, but never above the line. 

The blue, red and purple boxes contain combinations of 
{u, u) that can be reached for various choices of e, B, n, and 
Ub- The blue box covers the region in which B ^ lO'^^'^ G and 
n — 3, and the gravitational wave torque (Qi) dominates, 
i.e. />Gw 3> i'EM, where the subscripts EM and GW denote 
the electromagnetic and gravitational wave components of 
the spin down respectively. The red box covers the region 
in which the electromagnetic torque (Q2) dominates, with 
B = 10^"^ G and n — 3. The purple box also shows a region 
in which the Q2 term dominates, where we have chosen B = 
10^^ G and n — 2.3. As a rule of thumb, e determines the 
size of the box along the i^-axis, and i/t determines the size 
of the box along the />-axis. 

Let us first investigate what happens to the blue box 
when we vary the minimum and maximum ellipticity, emin 
and emax. The Qi term dominates in the region bounded by 
the blue box. The absolute value of the RQ slope increases 
as Emin dccrcascs, shrinking the range of z>. The curve PQ 
shifts to the left as emax increases, increasing |z>|, and hence 
lowering v. 

Let us now see what happens when we vary the min- 
imum and maximum magnetic field, Bmin and Bmax. The 
absolute value of the RS slope decreases as -Bmin in- 
creases, stretching the box sideways as we retreat from the 
gravitational- wave dominated limit. The blue box is always 
bounded above by the (n) = 5 age line. It shrinks, and flat- 
tens as the role of Qi diminishes. 

We now discuss the purple and red boxes in which Q2 



dominates. The region bounded by the purple box has B = 
10^"^ G, and n = 2.3, whereas the red box has the same 
B, but n = 3. Reducing n increases the spin-down rate by 
a factor of (7r7?*/c)"~^. Hence, for the same e and B, the 
purple box covers a smaller range of u than the red box. 
Both are considerably smaller than the blue box for the same 
range of e and Vb- Again, if emax increases, the purple and 
red boxes expand downwards. In Figure [2l we choose to plot 
the purple box with n — 2.3 because it lies partially within 
the sensitivity range of the search. Importantly, 1/ and u end 
up outside the search sensitivity range for n < 2.3 ox B > 
10^"^ G, restricting the range of astrophysical birth scenarios 
that our search is sensitive to. 

The range of v covered in the Q2-dominated limit is 
sensitive to B. In Figure[3l we show explicitly how varying B 
affects u, V. We plot eight red boxes, for 10^^ G (largest box) 
B ^ lO" '^ G (smallest box), and n = 3. As B mcreases, 
the red boxes shift to the left. For B ^ 5 x 10"G, the 
box falls out of the sensitivity range of the search. Also, the 
boxes shrink as Bmax increases. This happens because as B 
increases, i>EM increases. For B > 10^* G, we find v <^ Vb 
after 19 years, and the boxes end up on the (n) = 3 line. All 
the red boxes are bounded above by the (n) = 3 age line. 

Figures [5] and [3] provide constraints on the detectable 
range of e, B, n, over a broad range of Ub- We conclude that, 
in preparing to select the search templates, it is sensible 
to consider the parameter range lO"'^ ^ e, B ^ 10" G, 
2.3 ^ n ^ 5. A more detailed breakdown of the detectable 
and computationally feasible parameter ranges is presented 
in Section 16.21 Note that even though the particular boxes 
drawn as examples in Figures [2] and [3] do not cover the entire 
region between the sensitivity curve and the (n) = 1.8 line, 
one can potentially reach any point in that region with some 
combination of n, e and B. Also, each (:^, ;>) pair in the 
figures can be reached by an infinite set of combinations (e, 
B, n and Vb)- However, there are combinations of u and z> 
which are allowed in principle by age-based indirect limits 
but which cannot be reached from Vb with realistic choices 
of e, B, and n. 



6 TEMPLATE SPACING 

The cross-correlation search for SNR 1987A is computa- 
tionally limited rather than sensitivity limited over much 
of the parameter space in Figures [2] and [S] Therefore, the 
placement of templates is crucial. If the template grid is 
too coarse, the risk of missing the signal increases; if it 
is too fine, time is wasted searching redundant templates. 
In order to co mpute the optimal spacing, we construct 
a phase metric (|Balasubramanian et al.lll996l : lOwenI [l996l ) 
which computes the signal-to-noise ratio as a function of 
template spacing along each axis of the four-dimensional 
parameter space {u,e, B,n). The coherent phase metric for 
the conventional Taylor-expansion phase model is widely 
used in LIGO in bo th coherent and semi-coherent searches 
(Bra dv fc Creighton, 2000 ; Prix 2007 ; Wettc et al. 2008). al- 
though its semi-coherent form has not been fully investi- 
gated. In this section, we derive the semi-coherent phase 
metric for the astrophysical phase model defined by integrat- 
ing l|25p . We also estimate the range of detectable spin-down 
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v{Hz} 

Figure 2. Final states [u, i>) calculated from equation I I25II on the i^-li'l plane for a rang e of ellipticities (lO"'' ^ e ^ IQ-^), and birth 
spin frequencies (0.10 kHz i^t ^ 1.2 kHz), and for a 19 yr old pulsar. The blue lines surround the region where the Qi term dominates 
(B ^ lO^"*^'^ G, all n), the red lines surround the region where the Q2 term dominates {B = 10^^ G,n = 3), and the purple lines surround 
the region where the Q2 term dominates {B = 10^^ G, n = 2.3). The black curve shows the theoretical search sensitivity from solving 
equation (I17I I. The u <Si Uf, age limits are shown in cyan for (n) = 1.8 (top), (n) = 3 (middle) and (n) = 5 (bottom). 




v(Hz) 

Figure 3. Final states (f, i/) calculated from equation 12511 on the plane, for a range of magnetic field strengths. The eight red 

boxes surround regions which have n = 3 and cover the same range of e and vj, as Figure^ Their magnetic fields range from B = 10^^ G 
(largest box) to B = lO^*'^ G (smallest box). 
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values as well as magnetic field, ellipticity and braking index 
values given a computationally feasible number of templates. 

6.1 Semi-coherent phase metric 

When searching a template grid, it is extremely unlikely 
that one particular set of parameters will match the true 
signal exactly. What we have instead is a set of guessed 
parameters 6 + A0, describing the closest match, which are 
offset from the true values by a small amount, AO. For a 
given set of guessed parameters, the power spectrum of a 
time-coincident SFT pair is 



2A 



dt e 



iA#{t) 



(26) 



where A$(t) = ${t, 6 + AO) - ${t, 6) is the mismatch be- 
tween the actual and guessed phases, Tatart is the time at 
the midpoint of the SFT, and A is the gravitational wave 
amplitude. 

The mismatch between (I26|) and the power spectrum of 
the SFT pair if A0 = is defined to be 



m(0, A0) 



p{e,Ae) 



(27) 



P(0,O) 

and is related to the semi-coherent phase metric Sij by 

m{0,Ae) = s,j{0)Ae'Ae^, (28) 

where 1 ^ i,j ^ 4 label the various search parameters. 

For the cross-correlation search, we have — 
{i^,Qi,Q2,n). Hence, for a given mismatch m, the mini- 
mum (i.e. most conservative) template spacings are given 
by Au{0) = Vm/soo(6>) ,AQi(6>) = ^m/sii{e), AQ2{e) = 
\/rn/s22{0), An{Q) = \/ m/szz{Q). Note that it may be pos- 
sible to do better (i.e. expand the spacing) by taking advan- 
tage of the covariances between parameters embodied in the 
metric through (|28p : this issue deserves further study. 

In order to calculate Sij, we must first calculate the 
coherent phase metric, defined to be 



(29) 



with 



\ 1 fT'gtarl 



dt... and ftA$ = dA<!>/dA6' 
evaluated at AG = 0. Calculating gij analytically by inte- 
grating H25|) is non-trivial. However, a good approximation 
results if we integrate (|25p separately for the gravitational- 
wave and electromagnetic torques, and combine the answers 
in quadrature. Details of the calculation are shown in Ap- 
pendix A. In brief, tracking the gravitational-wave and elec- 
tromagnetic spin down separately yields two "sub-metrics", 
one comprising u and Qi (gravitational) and the other com- 
prising i/, Q2, and n (electromagnetic). Diagonal elements of 
Sij can be obtained by summing the two sub-metrics. 

The semi-coherent metric Sij is the average of the co- 
herent metric from Tstart = to Tstart = Tobs, where Tobs 
is the entire observation time spanned by all SFT pairs. It 
is defined to be Sij = {Tobs)~^ Jq°^^ dTsti^n gtj . From Ap- 
pendix B, the diagonal elements of the semi-coherent metric 
are: 



Soo 



_1_ 

36 



KiK2nu"+^QiQ2 



, 1 r^2 2 2n-2/ 

+ 360^^" " ^ 



(30) 



S22 



S33 



36 



'- lag-*^ obs 5 



36 



lag-^ obs 5 
2 ,.2nr 



36 



lag-' obs ■) 



(31) 
(32) 
(33) 



with Ki = Ki{v,Q\) and K2 = K2{i',Q2,ri). For pure 
gravitational-wave and electromagnetic spin down, we have 
{K\,K2) ~ (1,0) and (0,1) respectively. The full expres- 
sions for (|3U|l " (|33p are presented in Appendix B. Note that 
in p0|| - (|33|l . all frequency terms are normalised by v^ct- For 
clarity, we have set i/rcf ~ 1 Hz and do not display it. 

In Appendix C, we estimate the phase error which ac- 
cumulates after a time Tiag from mismatches in i^,Qi,Q2, 
and 71. We find that it scales with Tiag similarly to (|30p - 
(|33|l for Qi and Q2. For u, the phase error scales instead as 
Tiag, and for n, it scales as Tiag log [l -|- (n — l)Q2Tiagi^"~^] . 
In a semi-coherent search, the phase needs to be tracked to 
within 7r/4 over the interval Tag, not Tobs, unlike in fully 
coherent searches. Across the entire observation time Tobs, 
we require only that the frequency of the signal be tracked 
to within 1/AT. This adds an overall T^bs dependence to 



6.2 Computational cost of the search 

The run-time of the search code is proportional to 



where A'totai is the total number of templates 



required to search the parameter space. Trials with A^pairs = 
10^ comprising 1 year's worth of SFTs (from the two in- 
terferometers HI and LI), and Tag = Ihour take ~ Is per 
template on a single, 1-gigafiop computational node. We can 
therefore search ~ 10^ templates in a realistic run using 10^ 
nodes over two weeks. 

We now compare the computational cost of the as- 
trophysical phase model (|25|l against the Taylor-expansion 
model (|23|l . The semi-coherent metric for the latter model is 
not well studied, however recent work has yielded analytic 
expre ssions for the metric (Plctsch & Allen 2009; Plctsc^ 
I2OIOI ). Based on these expressions, we can estimate the num- 
ber of templates in the following way. 

Firstly, we consider the number of templates required 
to track the phase coherently over a time Tag. For the 
fc-th frequency derivative in the Taylor expansion model, 
the corresponding diagonal term of the c oherent metric 
scales as (Sij)^*^ oc T^f+^ (|Whitbeck|[200^ . The number 
of templates required to track the fc-th frequency deriva- 
tive coherently is then Nk oc 



»j;coh « The to- 

tal number of templates required for each coherent chunk 
of length Tag is therefore given by A^coh = 11^=0-'^'=' i-^- 
A^'coh oc ni—o-^iag^^i where s is the number of frequency 
derivatives required to track the gravitational wave phase 
(see Section[5ll. Now, assume that over a time Tobs, we sum a 
number of chunks incoherently, approximately proportional 
to Achunks oc Tobs /Tag B Now, usiug th e semi-coherent met- 
ric l|Pletsch fc Allenl2009l : |Pletschll2010l ). the number of tem- 



^ We emphasize that this is only an approximate estimate, as the 
cross-correlation method sums SFT pairs separated by a time up 
to and including Tiag. Strictly speaking, A^chunks > Tobs/Tiag. 
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plates required for s frequency derivatives is proportional to 
7s-^Vcoh, where 7^ is a 'refinement factor' which scales as 



chunks 

mately 



The total number of templates is then approxi- 



A^total OC N, 



s(3+l)/2 



/ j.^^ \ s(s + l)/2 



(34) 



(35) 



For the range of {vb, B,n) considered in Section l5.ll for 
Tlag ~ 1 hr, we must track terms up to and including z/'-'*^ 
in (|23p in order to keep the phase error overall below n/A. 
This gives iVtotai OC T^tTi^. 

Under the astrophysical phase model, we estimate 
A^totai = N^NqiNq2N„ from ((30l)-(l33|, where the sub- 
scripts denote the number of templates required for 
each individual parameter, e.g. N^, ~ v/Ai/. As (|30p 
yields different results for A^^ in the gravitational and 
electromagnetic limits, we bound A^,^ by taking it to 
be the sum of squares of the two limits, i.e. N^. = 

[Nl(Ki,K2)^(i.o)+Nl^Ki:K2)=(o.i)V^^- For a given mis- 
match m, we obtain 



Ar -2 10+3n 

iVtotal OC m ni' 



log{u)QiQ2 [Qi +n Q2 



^ lag-' obs 

(36) 

Equation (|36p is an approximate result, achieved by 
combining the two sub- metrics used in equations H30p - (|33p . 
It should be regarded as a rule of thumb. If gravitational- 
wave spin down dominates, we have S22 = S33 = 0, Tc = 
(4(51!/*)"^, and hence 



A''totai OC m V^r^ '^Ti T^^^. 



(37) 



If electromagnetic spin down dominates, we have sn — 0, 
Q2V^~^\ , and hence 



[(n- 



(38) 



The required template spacing therefore varies dramati- 
cally across the astrophysical parameter range. To illustrate, 
let us consider 0.1 kHz !^ IkHz, 10"^^ s^' < Qi ^ 



10 



s^ 10"^^ s^ ^ Q2 sC 10"" s^ and 2.3 ^ n sC 3.0, and 



hence 8 x 10"*^ ^ e ^ 8 x 10"*, 4xlO®G^B<4x lO" G. 
We assume a mismatch m of 0.2. The required resolutions 
in the four search parameters range across 

-4 



2.935 X 10"" s; Ai//Hz ^ 9.632 x 10' 

s; AQi/s^ ^ 3.685 x 10"^^ 
s; AQa/s^ ^ 2.120 x 10"^^ 



2.973 X 10" 
3.448 X 10" 
8.674 X 10 
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An 



^ 1.166 X 10* 



(39) 
(40) 
(41) 
(42) 



in this search volume. The number of templates required for 
each parameter is its range divided by its bin resolution. If 
the bin resolution is larger than its range, we require only 
one template. Equations (|39p - (|42|) imply a total number of 
templates between 2.958 x 10^ sC TVtotai 4.347 x 10^'' to 
cover the entire parameter space. Smaller values of v, Qi, Q2 
and n require fewer templates to cover their neighbourhood. 

Unfortunately, given the computational restrictions 
that we face, we cannot search the entire region of astro- 
physical parameters in Figure [21 In the following analysis, 
we therefore divide each axis in parameter space into (say) 
ten bins, i.e. a 10 x 10 x 10 x 10 hypercubic grid containing 



ut (kHz) 


e 


B(10" G) 


0.19-0.28 


> 1.6 X 10-4 


^ 2.0 


0.28-0.55 


> 1.0 X 10"'' 


$ 1-3 


0.55-1.00 


> 7.9 x 10-5 


$ 0.8 



Table 1. Table of 1^5, e, and B ranges (approximate) which are 
detectable by the cross-correlation search for SNR 1987A using 
LIGO S5 data, for 2.3 ^ n 3.0. The numbers in the table are 
based on the regions in which the computationally feasible (dark 
blue dots) and search-sensitive (cyan shaded) regions overlap in 
Figure O We assume standard values for the neutron star mass 
and radius, i.e. M* = 1.4M0, R^, = 10 km, and n = 3. 



10* "boxes", and calculate the localised resolution at the cen- 
troid of each box. The grid is spaced logarithmically along e 
and B to cover Qi and Q2 in a representative fashion. Only 
those boxes requiring A'^ < 10^ are practical to search. 

6.3 Astrophysical upper limits 

In this section, we combine the estimates of sensitivity and 
computational cost in Sectionsl4land [6^ respectivelv to iden- 
tify the ranges of the astrophysical parameters B and e that 
can be probed by a realistic search. In the event of a non- 
detection, upper limits on B and e can be placed. 

We solve (|25[l for a range of t'i,, e, and B, and calculate 
the characteristic wave strain ho from (|18[) . Figure |3] displays 
contours of ho versus B and e for n = 3 at two frequencies 
corresponding to = 300 Hz and vt = 1200 Hz. The cyan 
shaded areas indicate where ho ^ hth- The search is sensitive 
to a larger range of e and B as i/b rises. This occurs because 
the search sensitivity peaks at u ~ 150 Hz. For small vt and 
large e and B, the pulsar spins down after Tc = 19 yr to 
give v < 150 Hz. In the best case scenario, for = 1200 Hz, 
upper and lower limits on the magnetic field and ellipticity 
of B < 2.5 x 10" G and e > 8 X lO"'^ can be achieved. 

Unfortunately, the number of search templates required 
to cover the shaded region in Figure!?] is prohibitively large, 
as discussed in Section [6.21 Figure [5] shows both sensitivity 
and computational cost. Regions in which the search is sen- 
sitive (i.e. ho ^ hth) for given vi, and n are shaded in cyan. 
Overplotted as dark blue dots are the central coordinates 
of our grid boxes with A'^ ^ 10^. The panels correspond to 
a range of birth frequencies, Vb, and are grouped in pairs: 
n — 2.335 (left panel in pair) and n = 2.965 (right panel in 
pair) . The top pair shows the sensitivity and computational 
cost for Ub — 235 Hz, whereas the bottom pair corresponds 
to Ub = 955 Hz (bottom right). 

Figure [5] shows that the search sensitivity increases with 
i^b and e, and decreases with B. On the other hand, the com- 
putational efficiency of the search decreases with ft and e, 
and increases with B. Even so, there is substantial over- 
lap between the regions in which the search is sensitive and 
the regions which are computationally permissible. We note 
that as each individual dot in Figure [5] represents a region 
in which A^ ^ 10®, it is not feasible to search over all the 
dotted areas, as this would mean A^totai > 10^ Therefore, 
when implementing the search, we will choose an appropri- 
ate range of parameters such that A^totai ^ 10®, using Figure 
[5] as a guide. 

Table [1] summarises the approximate range of i/b, e, and 
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Figure 4. Contour plots of ho as a function of e and BilO^'^ G) for SNR 1987A for values of i^j = 250 (left panel) and 1200 Hz (right 
panel). We assume n = 3 and a pulsar age of 19 years, as the S5 run began in 2006. The cyan shaded areas correspond to ho ^ hjiii 
where h^^ is defined in II16I I. 
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Figure 5. Log-log contour plots of ho as a function of e and B (10*0) for birth frequencies Vf, = 235 Hz (top panel) and vi, = 955 Hz 
(bottom panel), and a range braking indices, n. The frequency of the signal u is obtained by solving II25I I. and integrating over = 19yr. 
The cyan shaded areas indicate the regions in which ho ^ ?tthi where hth is defined in I I16II . The panels are arranged in pairs. Each pair 
shows n = 2.335 (left) and n = 2.965 (right). The dark blue dots indicate parameter combinations for which one has A'^ ^ 10^. 
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B in which the two regions in Figure [5] overlap. If the pulsar 
in SNR 1987A was born with a frequency between 0.19 kHz 
and 0.28 kHz, the search is sensitive to e > 1.6 x 10~* and 
B < 2.0 X 10^^ G. This range narrows as vi, increases; for 
birth frequencies between 0.55 kHz and 1.00 kHz, the search 
is sensitive to e > 7.9 x lO"'' and B < 0.8 x lO" G. We note 
that these estimates, derived from the limits on Q\ and Q2, 
assume the standard values for the neutron star mass and ra- 
dius, Mi, = IAMq and i?* = 10 km. It is possible that SNR 
1987A contains a low-mass neutron star with ~ O.ISM© 
l|lmsheniklll"993 ). in which case the limit on the ellipticity, 
for Ub = 1.00 kHz, would be e > 2.5 x 10"". 

We now comment briefly on the relevance of these 
limits. The range of B listed in Table [T] is wit hin the 
expe c ted theoretical range d iscussed in Section [2] (|Michell 
ll994l : IOgelman fc Alpal 120041 '). The range of e listed in Ta- 
ble [T] however, is larger than the maximum ellipticity sus- 
tainable by the unmagnetized neutron star crust for many 
equations of state. For example, conventional neutron stars 
are expected to support e ^ 10~®, while hybrid quark- 
baryon or me son-condensate sta rs can support e ^ 10~^ 
(lUshomirskv et al. 2000; O wcni I2OO5I : iHorowitz fc Kadaul 
l200g[ ). However, some exotic models do allow for larger el- 
lipticities. Solid strange quark stars are predicted to be able 
to sustain e ^ 6 x 10' (|Owenll2005l). For low-mass neutron 
stars , the limit is e ^ 5 x 10" (|lmsheni^ll992l : iHorowitj 
I2OIOI I. We note also that these limits apply only to elasti- 
cally supported deform ations; rnagnetically supported defor- 
mations can be larg er (|Melatoj|2007l : lAkgiin fc WassermarJ 
I2OO8I : iHaskelllliooa l. Therefore, even placing the relatively 
large upper limit of e < 10~* on the putative neutron star 
in SNR 1987A will be useful to some degree in constraining 
its mass and/or equation of state. 



on the pulsar's ellipticity and magnetic field: e ^ 8 x 10~ , 
B ^ 2.0 X lO'^^ G. The search is also expected to be sensi- 
tive to electromagnetic braking indices 2.3 ^ n 5C 3.0. Its 
greatest weakness remains that it assumes n to be constant 
throughout the semi-coherent integration. Constant n is the 
simplest possible astrophysical scenario, and it certainly de- 
serves to be considered in its own right, in view of the over- 
whelming computational cost of a variable-n search. Never- 
theless, it is vital to recognize that the constant-n hypothesis 
covers a small fraction of the astrophysical parameter space. 

A search using gravitational wave data is anticipated to 
begin soon and would be the first application of the cross- 
correlation method to a continuous wave search. 
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7 CONCLUSION 

In this paper, we describe the steps taken to quantify the 
astrophysical significance of a cross-correlation search for the 
supernova remnant SNR 1987A in LIGO S5 data. 

• We estimate the theoretical sensitivity of the cross- 
correlation search, and compare it to the conservative, age- 
based, wave strain estimate. In the frequency band 75 Hz 
< < 450 Hz, the age-based estimate lies above the detec- 
tion threshold. 

• We introduce an alternative to the Taylor expansion 
model of the gravitational wave phase based on a set of four 
astrophysical search parameters {y,e, B,n). The new phase 
model renders a search for a neutron star like SNR 1987A 
with a high spin-down rate computationally feasible. 

• To estimate the optimal template spacing for the search, 
we calculate the semi-coherent phase metric corresponding 
to this astrophysical model. 

• We place detection limits on e and B for a range of birth 
spin frequencies, 0.1 kHz ^Ub^ 1.2 kHz. 

With the required template spacing and current com- 
putational capabilities discussed in Section 16.21 we will be 
able to search up to approximately 10^ templates. In the 
event of a non-detection, considering the parameter range 
discussed in this paper and assuming the standard neutron 
star mass and radius, we expect to place the following limits 
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APPENDIX A: CALCULATION OF THE COHERENT METRIC Gi 



This appendix details the calculation of the diagonal terms of the coherent metric, Qij (|29p . We start by evaluating the 
frequency v{t) at time t, by assuming that v(t) is a simple sum of separate, independent contributions from gravitational-wave 
and electromagnetic spin down: 



+ 



K2V 



(43) 
(44) 



(1 + AQivHf'^ [l + {n- l)Q2i^"-4]'''""' ■ 

Here, K\ and K2 are constants which satisfy K\ + K2 = 1, and the search parameters v,Q\,Q2, and n are defined at a 
reference time to. Recall that v is normalised by Urd, which we set to llfz and do not write down, for simplicity. The first 
term in (|44p follows from the first integral in (|43p by assuming Q2 = 0. The second term in (|44|) follows from the second 
integral in H43p by assuming Qi = 0. Needless to say, the exact solution for u(t) follows from solving H25p self-consistently 
for Qi 7^ 0, Q2 7^ 0, but this is too difficult to solve analytically. As the phase metric calculation is useful only in an analytic 
form, we adopt the approximation in (|43|l . 
The phase at time t is given by, 

nt + io 

dtu(t) 
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We can expand each term in the regimes {Qiv"^) ^ 3> t and {Q2V 
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Subtracting ^} from (|48)) gives 
A<l>(t) = Aiy{K2 + Ki) (^t + ^ - to 
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We now take the derivative of (|49p with respect to Av, AQ\, AQ2, and An. We have 
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We construct gij by substituting (|50p - (|53p into (I29p . In this paper, we only require the diagonal terms of the metric. The 
relevant terms (500,511,522,333) are 
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APPENDIX B: SEMI-COHERENT METRIC 

In this appendix, we list in full the diagonal terms of the semi-coherent metric presented in (|30p - (|33p . The relevant terms 

(soo, sii, S22, S33) are 
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APPENDIX C: ANALYTIC ACCURACY ESTIMATES FOR THE ASTROPHYSICAL PHASE MODEL 

In this appendix, we motivate (|30p - (|33p physically by calculating the phase error in two special cases: (i) pure gravitational- 
wave spin down, and (ii) pure electromagnetic spin down. In the gravitational wave case, (|25p reduces to 



dv „ 5 

u{t) = 



where we take i/rcf = 1 Hz for simplicity and = v{t = Q). The gravitational wave phase at t = Tiag is then 



$(riag) - $(to) = / dt,y{t) 

J to 

(1 + 4Ql/.*riag)^/* - 1 



(62) 
(63) 

(64) 
(65) 



There are two regimes to be considered: (i) Tiag ^ (4Qii/*)~^, and ii) Tiag <C (4Qi!/'*)~^. In terms of the characteristic 



4W{t)\- 



(66) 



the two regimes correspond to (i) Tiag ^ Tcito), and (ii) Tiag ^ Tc{to)- In the case of SNR 1987A, we have Tc ~ 19years (in 
2006, when the S5 run began) and Tag ~ Ihr, i.e. regime (ii). 
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Given small errors AQi and Ai/ in Qi and i/, the phase error that accumulates between the template and the signal after 
a time Tiag is 

A$ = — AO1 + — Ai/ 67 

dQ\ dv 

= -i/^'TiigAQi + TiagAi. (68) 

Overall, therefore, the number of templates required scales as Tj^g regardless of how rapidly the neutron star is spinning 
down. This scaling matches the conventional Taylor expansion if v and v suffice to track the signal (A'^totai oc Tjig) but is much 

more economical if v is needed (A^totai oc T[1^), which happens for u>l.7x 10-^ Hz (tc/IO^ yr) ^ (Tiag/l hour)"^ In the SNR 
1987A search, we cover frequencies above 0.1 kHz, so i> always contributes significantly. Hence the phase model 1)25^ is always 
preferable. 

Now suppose the electromagnetic term dominates. Equation (|25|1 with Urct ~ 1 Hz reduces to 

^ = ~Q2ry" (69) 

"(i) = TTT— TT (™) 

[l + (n-l)Q2i^"-4]'/^""" 

and the gravitational wave phase after a time Tiag is 

<E-(riag) = / dtu{t) (71) 

J to 

[1 + (n - 1) Qa^-^Tiag] - 1 
For small errors AQ2, Av and An in Q2, v and n, the phase error between the template at the signal after a time Tiag is 

Ci'I' A ^ A A /„oN 

A$ = — — AO2 + -;-Ai^ + — An (73) 
dQ2 dv dn 



'lag 



AQ2 + Ti^gAf + 



2 

Tiag log [1 + (n - 1) Q2riagi/""'] An (74) 

Hence in the electromagnetic limit, the phase error due to Af scales in the same way as in the gravitational wave limit. The 
phase error due to AQ2 scales as Tjig, and the phase error due to An scales as Tiag log [l + (n — 1) Q2Tlag'^"~^] • Overall, the 
number of templates required scales as Tjag log (Tiag). This represents a saving if the second frequency derivative is important 
which, as shown above, is true for the range of signal frequencies considered in this search. 
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